#include "atmospheric_correction.h"

// 蒙气改正
double MQC(double h) {
    return 0.0002967 / tan(h + 0.003138 / (h + 0.08919));
}

double MQC2(double ho) {
    return -0.0002909 / tan(ho + 0.002227 / (ho + 0.07679));
}

// 视差改正
void parallax(double z[3], double H, double fa, double high) {
    double dw = 1;
    if (z[2] < 500) dw = cs_AU; // 使用 astronomy_constants.h 中的 cs_AU
    z[2] *= dw;

    double r0, x0, y0, z0;
    double f = cs_ba; // 使用 astronomy_constants.h 中的 cs_ba
    double u = atan(f * tan(fa));
    double g = z[0] + H;

    r0 = cs_rEar * cos(u) + high * cos(fa); // 使用 astronomy_constants.h 中的 cs_rEar
    z0 = cs_rEar * sin(u) * f + high * sin(fa);
    x0 = r0 * cos(g);
    y0 = r0 * sin(g);

    double s[3];
    s[0] = z[0] * cos(z[1]) * cos(z[2]);
    s[1] = z[0] * cos(z[1]) * sin(z[2]);
    s[2] = z[0] * sin(z[1]);

    s[0] -= x0;
    s[1] -= y0;
    s[2] -= z0;

    z[0] = atan2(s[1], s[0]);
    z[1] = atan2(s[2], sqrt(s[0] * s[0] + s[1] * s[1]));
    z[2] = sqrt(s[0] * s[0] + s[1] * s[1] + s[2] * s[2]) / dw;
}